overall <- read_rds(here("data", "output_data", "overall.rds"))
black_miss <- read_rds(here("data", "output_data", "black.rds"))
black_complete <- read_rds(here("data", "output_data", "black.rds")) %>%
#removed states with partially missing data (i.e., between 4 and 19 yearly outcome data
#missing) #for four states (Alaska, New Mexico, Rhode Island, Utah) and completely
#missing #for nine states (Hawaii, Idaho, Maine, Montana, New Hampshire, North Dakota,
#South Dakota, Vermont, Wyoming. Total removed (n=4 + 9 = 13)
filter(!(state %in% c('Alaska', 'New Mexico','Rhode Island', 'Utah',
'Hawaii', 'South Dakota',
'Idaho', 'Maine', 'Montana', 'New Hampshire',
'North Dakota', 'Vermont', 'Wyoming')))
hispanic_miss <- read_rds(here("data", "output_data", "hispanic.rds"))
hispanic_complete <- read_rds(here("data", "output_data", "hispanic.rds")) %>%
#removed states with partially missing data (i.e., between 4 and 19 yearly outcome data missing)
# for twelve (12) states (Alabama, Arkansas, Idaho, Iowa, Kentucky, Louisiana,
# Minnesota, Nebraska, Rhode Island, South Carolina, Tennessee, Wisconsin) and completely missing for
# eleven (11) states (Alaska, Delaware, Maine, Mississippi, Montana, New Hampshire,
# North Dakota, South Dakota, Vermont, West Virginia, Wyoming). Total removed (n=11+12 = 23)
filter(!state %in% c("Alabama", "Arkansas", "Idaho", "Iowa", "Kentucky",
"Louisiana", "Minnesota", "Nebraska", "Rhode Island",
"South Carolina", "Tennessee", "Wisconsin",
"Alaska","Delaware","Maine","Mississippi","Montana",
"New Hampshire","North Dakota", "South Dakota",
"Vermont","West Virginia", "Wyoming")) %>%
#We imputed below outcome data for the six states with minimal missing (Kansas,
# Maryland, Missouri, North Carolina, Oregon, Utah). We imputed the
# outcome data with the closest outcome data (previous/next year) given
# that excluding these six states will have prevented the generalized
# synthetic control method from working
arrange(state, year) %>%
dplyr::group_by(state) %>%
fill(c("cvd_death_rate", "population"), .direction = "downup") %>%
dplyr::ungroup()
white <- read_rds(here("data", "output_data", "white.rds"))
men <- read_rds(here("data", "output_data", "men.rds"))
women <- read_rds(here("data", "output_data", "women.rds"))est_overall <- gsynth(cvd_death_rate ~ treatedpost + primarycare_rate +
cardio_rate + population +
low_educ + employed_for_wages + party +
low_income + married + male + race_nonwhite,
data = overall,
EM = F,
index = c("state_id","year"),
inference = "parametric", se = TRUE,
r = c(0, 5),
seed = 123,
nboots = 200, CV = TRUE, force = "two-way", parallel = FALSE)
#ATT plot
p1 <- plot(est_overall, type = "counterfactual", raw = "none",
theme.bw = TRUE, main = "", #ylim = c(120, 200),
legendOff = TRUE, xlab = "", ylab = "") +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
#axis.text.x = element_text(angle = 90, vjust = 0.5),
plot.margin = margin(10, 10, 50, 10))overall_by_period <- est_overall[["est.att"]] %>%
as.data.frame() %>%
rownames_to_column(var= "time") %>%
mutate(time=as.numeric(time)) %>%
dplyr::filter(time>=0)
overall_average <- data.frame(est_overall$est.avg)
overall_averageest_black <- gsynth(cvd_death_rate ~ treatedpost + primarycare_rate +
cardio_rate + population +
low_educ + employed_for_wages + party +
low_income + married + male,
data = black_complete,
EM = F,
index = c("state_id","year"),
inference = "parametric", se = TRUE,
r = c(0, 5),
seed = 123,
nboots = 200, CV = TRUE, force = "two-way", parallel = FALSE)
p2 <- plot(est_black, type = "counterfactual", raw = "none", theme.bw = TRUE, main = "",#ylim = c(220, 475),
legendOff = TRUE, xlab = "", ylab = "") +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
#axis.text.x = element_text(angle = 90, vjust = 0.5),
plot.margin = margin(10, 10, 50, 10))black_by_period <- est_black[["est.att"]] %>%
as.data.frame() %>%
rownames_to_column(var= "time") %>%
mutate(time=as.numeric(time)) %>%
dplyr::filter(time>=0)
black_average <- data.frame(est_black$est.avg)
black_averageest_hispanic <- gsynth(cvd_death_rate ~ treatedpost + primarycare_rate +
cardio_rate + population +
low_educ + employed_for_wages + party +
low_income + married + male,
data = hispanic_complete,
EM = F,
index = c("state_id","year"),
inference = "parametric", se = TRUE,
r = c(0, 5),
seed = 123,
nboots = 200, CV = TRUE, force = "two-way", parallel = FALSE)
p3 <- plot(est_hispanic, type = "counterfactual", raw = "none", theme.bw = TRUE, main = "",#ylim = c(220, 475),
legendOff = TRUE, xlab = "", ylab = "") +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
#axis.text.x = element_text(angle = 90, vjust = 0.5),
plot.margin = margin(10, 10, 50, 10))hispanic_by_period <- est_hispanic[["est.att"]] %>%
as.data.frame() %>%
rownames_to_column(var= "time") %>%
mutate(time=as.numeric(time)) %>%
dplyr::filter(time>=0)
hispanic_average <- data.frame(est_hispanic$est.avg)
hispanic_averageest_white <- gsynth(cvd_death_rate ~ treatedpost + primarycare_rate +
cardio_rate + population +
low_educ + employed_for_wages + party +
low_income + married + male,
data = white,
EM = F,
index = c("state_id","year"),
inference = "parametric", se = TRUE,
r = c(0, 5),
seed = 123,
nboots = 200, CV = TRUE, force = "two-way", parallel = FALSE)
p4 <- plot(est_white, type = "counterfactual", raw = "none", theme.bw = TRUE, main = "",#ylim = c(220, 475),
legendOff = TRUE, xlab = "", ylab = "") +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
#axis.text.x = element_text(angle = 90, vjust = 0.5),
plot.margin = margin(10, 10, 50, 10))white_by_period <- est_white[["est.att"]] %>%
as.data.frame() %>%
rownames_to_column(var= "time") %>%
mutate(time=as.numeric(time)) %>%
dplyr::filter(time>=0)
white_average <- data.frame(est_white$est.avg)
white_averageest_men <- gsynth(cvd_death_rate ~ treatedpost + primarycare_rate +
cardio_rate + population +
low_educ + employed_for_wages + party +
low_income + married + race_nonwhite,
data = men,
EM = F,
index = c("state_id","year"),
inference = "parametric", se = TRUE,
r = c(0, 5),
seed = 123,
nboots = 200, CV = TRUE, force = "two-way", parallel = FALSE)
p5 <- plot(est_men, type = "counterfactual", raw = "none", theme.bw = TRUE, main = "",#ylim = c(220, 475),
legendOff = TRUE, xlab = "", ylab = "") +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
#axis.text.x = element_text(angle = 90, vjust = 0.5),
plot.margin = margin(10, 10, 50, 10))men_by_period <- est_men[["est.att"]] %>%
as.data.frame() %>%
rownames_to_column(var= "time") %>%
mutate(time=as.numeric(time)) %>%
dplyr::filter(time>=0)
men_average <- data.frame(est_men$est.avg)
men_averageest_women <- gsynth(cvd_death_rate ~ treatedpost + primarycare_rate +
cardio_rate + population +
low_educ + employed_for_wages + party +
low_income + married + race_nonwhite,
data = women,
EM = F,
index = c("state_id","year"),
inference = "parametric", se = TRUE,
r = c(0, 5),
seed = 123,
nboots = 200, CV = TRUE, force = "two-way", parallel = FALSE)
p6 <- plot(est_women, type = "counterfactual", raw = "none", theme.bw = TRUE, main = "",#ylim = c(220, 475),
legendOff = TRUE, xlab = "", ylab = "") +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
#axis.text.x = element_text(angle = 90, vjust = 0.5),
plot.margin = margin(10, 10, 50, 10))women_by_period <- est_women[["est.att"]] %>%
as.data.frame() %>%
rownames_to_column(var= "time") %>%
mutate(time=as.numeric(time)) %>%
dplyr::filter(time>=0)
women_average <- data.frame(est_women$est.avg)
women_averageBaseline characheristics
## Table 1: Baseline Characteristics----
#overall <- read_rds(here("data", "output_data", "overall.rds"))
tab <- overall %>%
filter(year==2014) %>%
set_variable_labels(
cvd_death_rate = "CVD deaths per 100,000 persons aged 45-64 years",
primarycare_rate = "Primary care clinicians per 100,000 residents",
cardio_rate = "Cardiologists per 100,000 residents",
low_income = "Percentage of residents aged 45-64 years with annual income less than $15,000",
male = "Percentage of residents aged 45-64 years who are males",
race_nonwhite = "Percentage of residents aged 45-64 years who are Non-White",
married = "Percentage of residents aged 45-64 years who are married",
low_educ = "Percentage of residents aged 45-64 years without high school degree",
employed_for_wages = "Percentage of residents aged 45-64 years who are employed for wages",
party = "Percentage of political party") %>%
select(cvd_death_rate, primarycare_rate, cardio_rate,low_income,
male, race_nonwhite, married, low_educ, employed_for_wages,
party, treated) %>%
mutate_at(vars(low_educ, employed_for_wages,
low_income, married, male, race_nonwhite), multiply_by, e2=100) %>%
set_value_labels(treated = c("Medicaid Non-Expansion States"= 0,
"Medicaid Expansion States"=1),
party = c("Republican"=0, "Democrat"=1, "Split"=2)) %>%
modify_if(is.labelled, to_factor) %>%
tbl_summary(by = treated,
missing = 'no',
digits = list(all_continuous()~1,
all_categorical()~1),
statistic = list(all_continuous() ~ '{mean} ({sd})',
all_categorical() ~ '{n} ({p}%)')) %>%
modify_header(label = "**Characteristics**") %>%
add_overall() %>%
add_stat_label(location = "row") %>%
modify_spanning_header(c("stat_1","stat_2") ~ "**State Medicaid Expansion Status**") %>% bold_labels()
tab| Characteristics | Overall, N = 50 | State Medicaid Expansion Status | |
|---|---|---|---|
| Medicaid Non-Expansion States, N = 16 | Medicaid Expansion States, N = 34 | ||
| CVD deaths per 100,000 persons aged 45-64 years, Mean (SD) | 154.0 (47.0) | 177.0 (52.5) | 143.1 (40.5) |
| Primary care clinicians per 100,000 residents, Mean (SD) | 77.0 (12.8) | 67.6 (7.2) | 81.5 (12.5) |
| Cardiologists per 100,000 residents, Mean (SD) | 6.1 (2.2) | 5.2 (1.5) | 6.5 (2.3) |
| Percentage of residents aged 45-64 years with annual income less than $15,000, Mean (SD) | 9.6 (3.1) | 10.7 (3.3) | 9.1 (2.8) |
| Percentage of residents aged 45-64 years who are males, Mean (SD) | 42.4 (2.9) | 41.8 (3.1) | 42.7 (2.8) |
| Percentage of residents aged 45-64 years who are Non-White, Mean (SD) | 21.4 (12.9) | 23.2 (12.5) | 20.5 (13.3) |
| Percentage of residents aged 45-64 years who are married, Mean (SD) | 61.4 (4.4) | 61.9 (5.6) | 61.2 (3.8) |
| Percentage of residents aged 45-64 years without high school degree, Mean (SD) | 6.8 (3.1) | 8.0 (3.5) | 6.2 (2.8) |
| Percentage of residents aged 45-64 years who are employed for wages, Mean (SD) | 63.9 (6.9) | 61.8 (7.8) | 64.9 (6.2) |
| Percentage of political party, n (%) | |||
| Â Â Â Â Republican | 24.0 (48.0%) | 15.0 (93.8%) | 9.0 (26.5%) |
| Â Â Â Â Democrat | 13.0 (26.0%) | 0.0 (0.0%) | 13.0 (38.2%) |
| Â Â Â Â Split | 13.0 (26.0%) | 1.0 (6.2%) | 12.0 (35.3%) |
Overall and subgroup effect of the Medicaid expansion on CVD mortality
overall_average$group='overall'
black_average$group='black'
hispanic_average$group='hispanic'
white_average$group='white'
men_average$group='men'
women_average$group='women'
average_data <- rbind(overall_average, black_average, hispanic_average,
white_average, men_average, women_average)
average_data$group <- factor(average_data$group,
levels=c('overall', 'black', 'hispanic', 'white',
'men','women'))
average_table <- average_data %>%
transmute(Group = group,
`Adjusted Mean Difference (95%CI)`= paste(digits(Estimate,2),"", " " ,"(",
digits(CI.lower,2),", ",
digits(CI.upper,2), ")", sep=""))
knitr::kable(average_table, caption = "Overall and subgroup effect of the Medicaid expansion on CVD mortality") %>%
kableExtra::kable_styling(bootstrap_options = c("striped", "bordered"),
full_width = F)| Group | Adjusted Mean Difference (95%CI) | |
|---|---|---|
| ATT.avg | overall | -4.29 (-9.87, 1.29) |
| ATT.avg1 | black | -5.36 (-22.63, 11.91) |
| ATT.avg2 | hispanic | -4.28 (-30.08, 21.52) |
| ATT.avg3 | white | -3.18 (-8.30, 1.94) |
| ATT.avg4 | men | -5.96 (-15.42, 3.50) |
| ATT.avg5 | women | -3.34 (-8.05, 1.37) |
# tripple difference in difference
ddd_black_hispanic_vs_white_overall <- bind_rows(black_average,
hispanic_average,
white_average) %>%
dplyr::select(group, ATT=Estimate, CI.lower, CI.upper, se=S.E.) %>%
pivot_wider(., names_from=group,
values_from=c(ATT, CI.lower, CI.upper, se)) %>%
summarise(mean_black=round((ATT_black-ATT_white), 2),
se_black=sqrt(se_black^2 + se_white^2),
lower_black=round((mean_black-qnorm(p=1-(0.05)/2, mean=0, sd=1)*se_black), 2),
upper_black=round((mean_black+qnorm(p=1-(0.05)/2, mean=0, sd=1)*se_black), 2),
mean_hispanic=round((ATT_hispanic-ATT_white), 2),
se_hispanic=sqrt(se_hispanic^2 + se_white^2),
lower_hispanic=round((mean_hispanic-qnorm(p=1-(0.05)/2, mean=0, sd=1)*se_hispanic), 2),
upper_hispanic=round((mean_hispanic+qnorm(p=1-(0.05)/2, mean=0, sd=1)*se_hispanic), 2)) %>%
pivot_longer(cols=everything(),
names_to = c(".value", "race"),
names_sep="_") %>%
arrange(race)
ddd_women_vs_men_overall <- bind_rows(men_average,
women_average) %>%
dplyr::select(group, ATT=Estimate, CI.lower, CI.upper, se=S.E.) %>%
pivot_wider(., names_from=group,
values_from=c(ATT, CI.lower, CI.upper, se)) %>%
summarise(mean=round((ATT_women-ATT_men), 2),
se=sqrt(se_men^2 + se_women^2),
lower=round((mean-qnorm(p=1-(0.05)/2, mean=0, sd=1)*se), 2),
upper=round((mean+qnorm(p=1-(0.05)/2, mean=0, sd=1)*se), 2))
bind_rows(ddd_black_hispanic_vs_white_overall, ddd_women_vs_men_overall) %>%
rename(group = race) %>%
mutate(group = case_when(group=="black"~"Black vs White",
group=="hispanic"~"Hispanic vs White",
TRUE~"Women vs men")) %>%
mutate(`Adjusted Difference in Mean Difference (95%CI)`= paste(digits(mean,2),"", " " ,"(",
digits(lower,2),", ",
digits(upper,2), ")", sep="")) %>%
dplyr::select(Group= group,`Adjusted Difference in Mean Difference (95%CI)`) %>%
knitr::kable(caption = "Difference in mean differences") %>%
kableExtra::kable_styling(bootstrap_options = c("striped", "bordered"),
full_width = F)| Group | Adjusted Difference in Mean Difference (95%CI) |
|---|---|
| Black vs White | -2.18 (-20.19, 15.83) |
| Hispanic vs White | -1.10 (-27.40, 25.20) |
| Women vs men | 2.62 (-7.95, 13.19) |
Analytical data structure of the Medicaid expansion states and control states for the overall study population.
q1 <- panelview(cvd_death_rate ~treatedpost , data=overall, index=c("state","year"),
pre.post= TRUE,
cex.axis.x=20,
cex.axis.y=20,
cex.lab=20)Overall annual effect of Medicaid expansion on CVD deaths per 100,000 persons overall and by race and sex.
f_overall <- overall_by_period %>%
ggplot(aes(x = time, y = ATT)) +
geom_point(color="#030303") +
geom_errorbar(aes( ymax = CI.lower, ymin = CI.upper) , width = .2, size = 0.7,
color="#030303")+
geom_line(color="#030303") +
ylim(-38, 38)+
scale_x_continuous(breaks=seq(0,6,1))+
geom_hline(yintercept = 0, lty=2,color="#030303") +
ggtitle('Overall')+
labs(x = "Number of years since Medicaid Expansion",
y = "Change in CVD mortality, per 100,000 population") +
theme_classic(base_size = 15) +
theme(axis.line = element_line(colour = "grey50", size = 1),
panel.background = element_blank(),
axis.title = element_blank(),
plot.title = element_text(hjust = 0.5))
f_black_hispanic_white <- bind_rows(black_by_period %>% mutate(group="black"),
hispanic_by_period %>% mutate(group="hispanic"),
white_by_period %>% mutate(group="white")) %>%
ggplot(aes(x = time, y = ATT, color = group)) +
scale_color_manual(values=c("#8A2BE2", "#FF7F50", "#458B00"),
labels=c('Black', 'Hispanic','White'))+
geom_point(position = position_dodge(width = .5)) +
geom_errorbar( aes( ymax = CI.upper, ymin = CI.lower) , width = .2,
position = position_dodge(width = .5), size = 0.7 )+
geom_line(position = position_dodge(width = .5), size = 0.7) +
ylim(-38, 38)+
scale_x_continuous(breaks=seq(0,6,1))+
#ggtitle('B')+
geom_hline(yintercept = 0, lty=2, color="black") +
theme(legend.position = "top",
axis.title = element_blank(),
axis.line = element_line(colour = "grey50", size = 1),
panel.background = element_blank(),
legend.background = element_rect(fill = "lemonchiffon",
colour = "grey50",
size = 1),
legend.title = element_blank())
f_men_women <- bind_rows(men_by_period %>% mutate(group="men"),
women_by_period %>% mutate(group="women")) %>%
ggplot(aes(x = time, y = ATT, color = group)) +
scale_color_manual(values=c("#1874CD", "#EE3B3B"),
labels=c('Men','Women'))+
geom_point(position = position_dodge(width = .5)) +
geom_errorbar( aes( ymax = CI.upper, ymin = CI.lower) , width = .2,
position = position_dodge(width = .5), size = 0.7 )+
geom_line(position = position_dodge(width = .5), size = 0.7) +
ylim(-38, 38)+
scale_x_continuous(breaks=seq(0,6,1))+
#ggtitle('B')+
geom_hline(yintercept = 0, lty=2, color="black") +
theme(legend.position = "top",
axis.title = element_blank(),
axis.line = element_line(colour = "grey50", size = 1),
panel.background = element_blank(),
legend.background = element_rect(fill = "lemonchiffon",
colour = "grey50",
size = 1),
legend.title = element_blank())
f_by_overall_race_sex <- ggarrange(f_overall, f_black_hispanic_white, f_men_women,
ncol = 1, nrow=3) %>%
annotate_figure(left = textGrob("Change in CVD mortality, per 100,000 population",
rot = 90, vjust = 1, gp = gpar(cex = 1.3)), #cex = 0.8
bottom = textGrob("Number of years since Medicaid Expansion", gp = gpar(cex = 1.3)))
f_by_overall_race_sex
## Figure 3. Pre-treatment fit
Pre-treatment fit showing the observed CVD deaths per 100,000
persons along with the counterfactual (synthetic control). Vertical line
represents the beginning of expansion. The solid line is expansion
state, dashed line is the synthetic control.
p1_6_fit <- ggarrange(
p1+ rremove("ylab") + rremove("xlab")+ggtitle("Overall"),
p2+ rremove("ylab") + rremove("xlab")+ggtitle("Black"),
p3+ rremove("ylab") + rremove("xlab")+ggtitle("Hispanic"),
p4+ rremove("ylab") + rremove("xlab")+ggtitle("White"),
p5+ rremove("ylab") + rremove("xlab")+ggtitle("Men"),
p6+ rremove("ylab") + rremove("xlab")+ggtitle("Women"),
font.label = list(size = 2, color = "black", face = "bold", family = NULL),
ncol=2, nrow=3,common.legend = TRUE)
p1_6_fitOverall annual effect of the Medicaid expansion on CVD deaths per 100,000 persons for overall population and for the Black, Hispanic, White, Men and women subpopulations.
average_by_period_table <-
bind_rows(overall_by_period %>% mutate(group="Overall"),
black_by_period %>% mutate(group="Black"),
hispanic_by_period %>% mutate(group="Hispanic"),
white_by_period %>% mutate(group="White"),
men_by_period %>% mutate(group="Men"),
women_by_period %>% mutate(group="Women")) %>%
mutate(`Adjusted Mean Difference (95%CI)`= paste(digits(ATT,2),"", " " ,"(",
digits(CI.lower,2),", ",
digits(CI.upper,2), ")", sep="")) %>%
dplyr::select(group, time, `Adjusted Mean Difference (95%CI)`)
average_by_period_table %>%
knitr::kable(caption = "Overall annual effect of the Medicaid expansion on CVD deaths per 100,000 persons for the overall population and for the Black, Hispanic, White, Men and women subpopulations") %>%
kableExtra::kable_styling(bootstrap_options = c("striped", "bordered"),
full_width = F)| group | time | Adjusted Mean Difference (95%CI) |
|---|---|---|
| Overall | 0 | 0.40 (-3.38, 4.17) |
| Overall | 1 | -5.09 (-9.34, -0.85) |
| Overall | 2 | -7.61 (-13.43, -1.79) |
| Overall | 3 | -2.81 (-9.33, 3.71) |
| Overall | 4 | -4.18 (-10.76, 2.41) |
| Overall | 5 | -3.98 (-11.23, 3.28) |
| Overall | 6 | -1.55 (-9.88, 6.77) |
| Black | 0 | -1.48 (-16.15, 13.19) |
| Black | 1 | -6.19 (-23.10, 10.72) |
| Black | 2 | -14.67 (-32.74, 3.40) |
| Black | 3 | -4.66 (-23.60, 14.29) |
| Black | 4 | -5.44 (-26.24, 15.37) |
| Black | 5 | -1.61 (-25.21, 22.00) |
| Black | 6 | 1.45 (-24.93, 27.82) |
| Hispanic | 0 | -2.56 (-18.05, 12.93) |
| Hispanic | 1 | -6.31 (-27.21, 14.60) |
| Hispanic | 2 | -7.93 (-31.09, 15.23) |
| Hispanic | 3 | -4.26 (-31.20, 22.68) |
| Hispanic | 4 | 0.32 (-28.57, 29.21) |
| Hispanic | 5 | -6.78 (-37.26, 23.70) |
| Hispanic | 6 | -0.18 (-34.13, 33.78) |
| White | 0 | 0.02 (-3.63, 3.68) |
| White | 1 | -5.42 (-10.99, 0.15) |
| White | 2 | -4.48 (-10.09, 1.13) |
| White | 3 | -2.16 (-8.10, 3.77) |
| White | 4 | -1.88 (-7.63, 3.87) |
| White | 5 | -3.06 (-9.81, 3.69) |
| White | 6 | -1.67 (-9.49, 6.15) |
| Men | 0 | 1.50 (-5.36, 8.36) |
| Men | 1 | -4.72 (-12.59, 3.15) |
| Men | 2 | -10.35 (-20.28, -0.42) |
| Men | 3 | -2.40 (-12.92, 8.13) |
| Men | 4 | -8.16 (-19.81, 3.49) |
| Men | 5 | -5.74 (-18.13, 6.65) |
| Men | 6 | -4.16 (-17.59, 9.28) |
| Women | 0 | -1.58 (-6.37, 3.20) |
| Women | 1 | -5.50 (-9.92, -1.09) |
| Women | 2 | -4.36 (-9.87, 1.15) |
| Women | 3 | -2.91 (-8.15, 2.34) |
| Women | 4 | -2.06 (-7.69, 3.58) |
| Women | 5 | -4.57 (-10.57, 1.44) |
| Women | 6 | -0.08 (-5.49, 5.34) |
Sensitivity analysis evaluating the effect of using heterogenous samples vs homogenous samples on both adjusted mean differences and difference in mean differences.
#Whites with black sample (n=37 states included)
white_withblacksample <- read_rds(here("data", "output_data", "white.rds")) %>%
filter(!(state %in% c('Alaska', 'New Mexico','Rhode Island', 'Utah',
'Hawaii', 'South Dakota',
'Idaho', 'Maine', 'Montana', 'New Hampshire',
'North Dakota', 'Vermont', 'Wyoming')))
est_white_withblacksample <- gsynth(cvd_death_rate ~ treatedpost + primarycare_rate +
cardio_rate + population +
low_educ + employed_for_wages + party +
low_income + married + male,
data = white_withblacksample,
EM = F,
index = c("state_id","year"),
inference = "parametric", se = TRUE,
#min.T0 = 5,
r = c(0, 5),
seed = 123,
nboots = 200, CV = TRUE, force = "two-way", parallel = FALSE)## Cross-validating ...
## r = 0; sigma2 = 41.47336; IC = 4.63849; PC = 34.56113; MSPE = 52.04339
## r = 1; sigma2 = 28.19493; IC = 4.91483; PC = 133.70899; MSPE = 37.70570*
## r = 2; sigma2 = 21.85741; IC = 5.27680; PC = 189.27671; MSPE = 40.22883
## r = 3; sigma2 = 17.07014; IC = 5.60049; PC = 214.83195; MSPE = 43.51960
## r = 4; sigma2 = 15.45230; IC = 6.02614; PC = 255.26004; MSPE = 49.83340
## r = 5; sigma2 = 13.26928; IC = 6.35339; PC = 271.50968; MSPE = 74.10160
##
## r* = 1
##
##
Simulating errors .....
Bootstrapping ...
## ..
est_white_withblacksample_average <- data.frame(est_white_withblacksample$est.avg)
est_white_withblacksample_average %>%
knitr::kable(caption = "Sensitivity analysis for the adjusted mean difference: Whites and Blacks have same sample (n=37)") %>%
kableExtra::kable_styling(bootstrap_options = c("striped", "bordered"),
full_width = F)| Estimate | S.E. | CI.lower | CI.upper | p.value | |
|---|---|---|---|---|---|
| ATT.avg | -2.15793 | 3.740264 | -9.488713 | 5.172853 | 0.563976 |
#Whites with Hispanic sample (n=27 states included)
white_withhispanicsample <- read_rds(here("data", "output_data", "white.rds")) %>%
filter(!state %in% c("Alabama", "Arkansas", "Idaho", "Iowa", "Kentucky",
"Louisiana", "Minnesota", "Nebraska", "Rhode Island",
"South Carolina", "Tennessee", "Wisconsin",
"Alaska","Delaware","Maine","Mississippi","Montana",
"New Hampshire","North Dakota", "South Dakota",
"Vermont","West Virginia", "Wyoming"))
est_white_withhispanicsample<- gsynth(cvd_death_rate ~ treatedpost + primarycare_rate +
cardio_rate + population +
low_educ + employed_for_wages + party +
low_income + married + male,
data = white_withhispanicsample,
EM = F,
index = c("state_id","year"),
inference = "parametric", se = TRUE,
#min.T0 = 5,
r = c(0, 5),
seed = 123,
nboots = 200, CV = TRUE, force = "two-way", parallel = FALSE)## Cross-validating ...
## r = 0; sigma2 = 27.87331; IC = 4.56308; PC = 20.90499; MSPE = 55.58006*
## r = 1; sigma2 = 16.13953; IC = 4.86382; PC = 141.06632; MSPE = 60.82607
## r = 2; sigma2 = 12.59613; IC = 5.39248; PC = 210.92380; MSPE = 63.71601
## r = 3; sigma2 = 7.25834; IC = 5.54719; PC = 179.74642; MSPE = 73.55888
## r = 4; sigma2 = 4.99104; IC = 5.80804; PC = 163.69302; MSPE = 76.60868
## r = 5; sigma2 = 2.24814; IC = 5.57526; PC = 91.82514; MSPE = 122.77355
##
## r* = 0
##
##
Simulating errors .....
Bootstrapping ...
## ..
est_white_withhispanicsample_average <- data.frame(est_white_withhispanicsample$est.avg)
est_white_withhispanicsample_average %>%
knitr::kable(caption = "Sensitivity analysis for the adjusted mean difference: Whites and Hispanics have same sample (n=27)") %>%
kableExtra::kable_styling(bootstrap_options = c("striped", "bordered"),
full_width = F)| Estimate | S.E. | CI.lower | CI.upper | p.value | |
|---|---|---|---|---|---|
| ATT.avg | -11.44491 | 3.651427 | -18.60158 | -4.288245 | 0.0017223 |
ddd_black_vs_white_overall_homosample <- bind_rows(black_average %>%
mutate(group="black"),
est_white_withblacksample_average %>%
mutate(group="white")) %>%
dplyr::select(group, ATT=Estimate, CI.lower, CI.upper, se=S.E.) %>%
pivot_wider(., names_from=group,
values_from=c(ATT, CI.lower, CI.upper, se)) %>%
summarise(mean=round((ATT_black-ATT_white), 2),
se=sqrt(se_black^2 + se_white^2),
lower=round((mean-qnorm(p=1-(0.05)/2, mean=0, sd=1)*se), 2),
upper=round((mean+qnorm(p=1-(0.05)/2, mean=0, sd=1)*se), 2),
pvalue=round(pnorm(0,mean=abs(mean),sd=se) +
(1 - pnorm(0,mean=-abs(mean),sd=se)),2))
ddd_black_vs_white_overall_homosample %>%
knitr::kable(caption = "Sensitivity analysis for the difference in mean differences: Whites and Blacks have same sample (n=37)") %>%
kableExtra::kable_styling(bootstrap_options = c("striped", "bordered"),
full_width = F)| mean | se | lower | upper | pvalue |
|---|---|---|---|---|
| -3.2 | 9.570525 | -21.96 | 15.56 | 0.74 |
ddd_hispanic_vs_white_overall_homosample <- bind_rows(hispanic_average %>%
mutate(group="hispanic"),
est_white_withhispanicsample_average %>% mutate(group="white")) %>%
dplyr::select(group, ATT=Estimate, CI.lower, CI.upper, se=S.E.) %>%
pivot_wider(., names_from=group,
values_from=c(ATT, CI.lower, CI.upper, se)) %>%
summarise(mean=round((ATT_hispanic-ATT_white), 2),
se=sqrt(se_hispanic^2 + se_white^2),
lower=round((mean-qnorm(p=1-(0.05)/2, mean=0, sd=1)*se), 2),
upper=round((mean+qnorm(p=1-(0.05)/2, mean=0, sd=1)*se), 2),
pvalue=round(pnorm(0,mean=abs(mean),sd=se) +
(1 - pnorm(0,mean=-abs(mean),sd=se)),2))
ddd_hispanic_vs_white_overall_homosample %>%
knitr::kable(caption = "Sensitivity analysis for the difference in mean differences: Whites and Hispanics have same sample (n=27)") %>%
kableExtra::kable_styling(bootstrap_options = c("striped", "bordered"),
full_width = F)| mean | se | lower | upper | pvalue |
|---|---|---|---|---|
| 7.17 | 13.65951 | -19.6 | 33.94 | 0.6 |
Analytical data structure of Medicaid expansion states and control states for Blacks (before and after removing states with missing outcomes.
qblack_miss <- panelview(cvd_death_rate ~treatedpost , data=black_miss, index=c("state","year"), pre.post= TRUE)qblack_complete <- panelview(cvd_death_rate ~treatedpost , data=black_complete, index=c("state","year"), pre.post= TRUE)ggarrange(qblack_miss+ rremove("ylab") + rremove("xlab")+ggtitle("Black (with missing data)"),
qblack_complete+ rremove("ylab") + rremove("xlab")+ggtitle("Black (with complete data)"),
ncol=1, nrow=2,
common.legend = TRUE)Analytical data structure of Medicaid expansion states and control states for Hispanics (before and after removing states with missing outcomes and nearest year imputation).
qhispanic_miss <- panelview(cvd_death_rate ~treatedpost ,
data=hispanic_miss, index=c("state","year"),
pre.post= TRUE)qhispanic_complete <- panelview(cvd_death_rate ~treatedpost , data=hispanic_complete, index=c("state","year"), pre.post= TRUE)ggarrange(qhispanic_miss+ rremove("ylab") + rremove("xlab")+ggtitle("Hispanic (with missing data)"),
qhispanic_complete+ rremove("ylab") + rremove("xlab")+ggtitle("Hispanic (with complete and imputed data)"),
ncol=1, nrow=2,
common.legend = TRUE)Annual difference in mean difference between the effect of the Medicaid expansion on CVD deaths per 100,000 persons
## Black or Hispanic vs White-
ddd_black_hispanic_vs_white_by_period <- bind_rows(black_by_period %>% mutate(group="black"),
hispanic_by_period %>% mutate(group="hispanic"),
white_by_period %>% mutate(group="white")) %>%
dplyr::select(time, group, ATT, CI.lower, CI.upper, se=S.E.) %>%
pivot_wider(., names_from=group,
values_from=c(ATT, CI.lower, CI.upper, se)) %>%
summarise(mean_black=round((ATT_black-ATT_white), 2),
se_black=sqrt(se_black^2 + se_white^2),
lower_black=round((mean_black-qnorm(p=1-(0.05)/2, mean=0, sd=1)*se_black), 2),
upper_black=round((mean_black+qnorm(p=1-(0.05)/2, mean=0, sd=1)*se_black), 2),
mean_hispanic=round((ATT_hispanic-ATT_white), 2),
se_hispanic=sqrt(se_hispanic^2 + se_white^2),
lower_hispanic=round((mean_hispanic-qnorm(p=1-(0.05)/2, mean=0, sd=1)*se_hispanic), 2),
upper_hispanic=round((mean_hispanic+qnorm(p=1-(0.05)/2, mean=0, sd=1)*se_hispanic), 2),
time=time) %>%
relocate(time) %>% pivot_longer(cols=-time,
names_to = c(".value", "race"),
names_sep="_") %>%
arrange(race)
s1 <- ddd_black_hispanic_vs_white_by_period %>%
ggplot(aes(x = time, y = mean, color = race)) +
scale_color_manual(values=c('#e41a1c', '#377eb8'), labels=c('Black vs White', 'Hispanic vs White'))+
geom_point(position = position_dodge(width = .5)) +
geom_errorbar( aes( ymax = upper, ymin = lower) , width = .2, position = position_dodge(width = .5), size = 0.7 )+
geom_line(position = position_dodge(width = .5), size = 0.7) +
ylim(-38, 38)+
scale_x_continuous(breaks=seq(0,6,1))+
geom_hline(yintercept = 0, lty=2) +
theme(legend.position = "top",
axis.title = element_blank(),
axis.line = element_line(colour = "grey50", size = 1),
panel.background = element_blank(),
legend.background = element_rect(fill = "lemonchiffon",
colour = "grey50",
size = 1),
legend.title = element_blank())
## Women vs men
ddd_women_vs_men_by_period <- bind_rows(men_by_period %>% mutate(group="men"),
women_by_period %>% mutate(group="women")) %>%
dplyr::select(time, group, ATT, CI.lower, CI.upper, se=S.E.) %>%
pivot_wider(., names_from=group,
values_from=c(ATT, CI.lower, CI.upper, se)) %>%
summarise(mean=round((ATT_women-ATT_men), 2),
se=sqrt(se_men^2 + se_women^2),
lower=round((mean-qnorm(p=1-(0.05)/2, mean=0, sd=1)*se), 2),
upper=round((mean+qnorm(p=1-(0.05)/2, mean=0, sd=1)*se), 2),
pvalue=round(pnorm(0,mean=abs(mean),sd=se) +
(1 - pnorm(0,mean=-abs(mean),sd=se)),2),
time=time) %>%
relocate(time)
s2 <- ddd_women_vs_men_by_period %>%
ggplot(aes(x = time, y = mean, color = '#4daf4a')) +
scale_color_manual(values=c('#4daf4a'), labels=c('Women vs Men'))+
geom_point(position = position_dodge(width = .5)) +
geom_errorbar( aes( ymax = upper, ymin = lower) , width = .2, position = position_dodge(width = .5), size = 0.7 )+
geom_line(position = position_dodge(width = .5), size = 0.7) +
ylim(-38, 38)+
geom_hline(yintercept = 0, lty=2) +
scale_x_continuous(breaks=seq(0,6,1))+
theme(legend.position = "top",
axis.line = element_line(colour = "grey50", size = 1),
panel.background = element_blank(),
legend.background = element_rect(fill = "lemonchiffon",
colour = "grey50",
size = 1),
legend.title = element_blank(),
axis.title = element_blank())
triple_ddd <- ggarrange(s1,s2, ncol = 2, nrow=1) %>%
annotate_figure(left = textGrob("Change in CVD mortality, per 100,000 population",
rot = 90, vjust = 1, gp = gpar(cex = 1.3)), #cex = 0.8
bottom = textGrob("Number of years since Medicaid Expansion", gp = gpar(cex = 1.3)))
triple_ddd